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SUMMARY 


For  coupled  structural  acoustic  scattering  problems,  the  spectral  response  often 
varies  rapidly  over  a  large  dynamic  range,  requiring  the  response  to  be  calculated  at 
many  closely  spaced  frequencies.  In  such  problems,  the  computation  time  is  often 
dominated  by  the  calculation  of  acoustic  interactions  which  do  not  vary  rapidly  with 
frequency.  It  is  shown  in  this  report  that  large  savings  in  computation  time  can 
be  achieved  without  significant  loss  of  accuracy  by  interpolating  on  appropriately 
scaled  intermediate  acoustic  quantities.  These  quantities  are  computed  on  a  coarse 
frequency  grid  and  then  interpolated  to  any  desired  frequency  resolution  in  order  to 
obtain  the  final  response.  The  overall  savings  in  computation  time  is  typically  more 
than  an  order  of  magnitude. 
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INTRODUCTION 


A  technique  for  solving  acoustic-structure  interaction  problems  was  described 
in  a  previous  Naval  Ocean  Systems  Center  (NOSC)  technical  report  (Schenck  and 
Benthien,  1989).  In  this  technique  the  elastic  structure  is  modeled  using  the  finite 
element  method,  and  the  acoustic  field  is  modeled  using  the  Combined  Helmholtz 
Integral  Equation  Formulation  (CHIEF)  (Schenck,  1968).  The  purpose  of  this  re¬ 
port  is  to  describe  a  frequency  interpolation  scheme  for  improving  the  computational 
efficiency  of  this  numerical  technique. 

Due  to  the  resonant  nature  of  the  elastic  structure,  it  is  often  necessary  to  compute 
the  acoustic  field  at  a  large  number  of  closely  spaced  frequencies.  In  such  problems 
the  computational  time  is  dominated  by  the  acoustic  calculations.  It  will  be  shown 
that  the  computational  time  can  be  reduced  considerably  by  interpolating  on  certain 
scaled  acoustic  quantities  which  vary  slowly  with  frequency.  Furthermore,  it  will  be 
shown  that  this  frequency  interpolation  scheme  results  in  little,  if  any,  loss  of  accuracy. 

The  second  section  of  this  report  contains  a  brief  summary  of  the  basic  equations 
involved  in  the  numerical  solution  of  acoustic-structure  interaction  problems.  The 
third  section  contains  a  description  of  the  frequency  interpolation  scheme  as  it  is  ap¬ 
plied  to  certain  appropriately  scaled  acoustic  quantities.  The  fourth  section  contains 
results  from  two  example  problems  which  illustrate  the  accuracy  and  computational 
efficiency  of  the  frequency  interpolation  scheme.  Such  improvements  in  computa¬ 
tional  efficiency  should  be  very  important  as  one  attempts  to  employ  this  numerical 
technique  at  higher  frequencies  where  computation  times  can  be  quite  large.  The  last 
section  summarizes  the  results  of  this  work. 
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BASIC  EQUATIONS 


Throughout  this  section  all  pressures,  displacements,  and  velocities  will  be  repre¬ 
sented  by  their  complex,  frequency  dependent  Fourier  components  corresponding  to 
the  time  dependence  e,ult.  The  elastic  structure  is  modeled  using  the  finite  element 
equations  (Zienkiewicz,  1977) 


(-u2M  +  K)U  =  FTad  +  Fd,  (1) 

where  M  is  the  mass  matrix,  K  is  the  stiffness  matrix,  FTad  is  the  load  vector  cor¬ 
responding  to  acoustic  loading,  FD  is  the  load  vector  corresponding  to  nonacoustic 
mechanical  drive  forces  on  the  structure,  and  U  is  a  column  vector  whose  components 
are  the  displacement  degrees-of- freedom  at  all  the  nodes  in  the  body.  The  components 
of  FTad  are  given  by 

K;'  =  -jsP<t>mndS,  (2) 

where  p  is  the  acoustic  pressure,  n  is  a  unit  normal  to  the  wet-surface  S  pointing  into 
the  fluid,  and  (f>m  is  a  finite-element  vector  interpolation  function.  The  displacement 
u(x)  at  a  point  x  in  the  body  can  be  expressed  in  terms  of  the  interpolation  functions 
4>m  as  follows: 

M 

u(x)  =  £  Um<t>m(x),  (3) 

771=1 

where  ... ,  Um  are  the  components  of  U. 

The  CHIEF  formulation  for  the  acoustic  field  is  based  on  the  Helmholtz  integral 
relations.  The  scalar  acoustic  wave  equation  for  the  acoustic  medium  is  replaced  by 
an  integral  equation  over  the  wet-surface  S  which  in  turn  is  approximated  by  a  system 
of  algebraic  equations  (Benthien,  Barach,  and  Gillette,  1988) 

AP  =  BV  +  Pinc ,  (4) 

where  P  and  V  are  column  vectors  whose  nth  components  are  the  pressure  and  normal 
velocity  (assumed  to  be  constant)  on  the  nth  subdivision  Sn  of  S’;  and  P,nc  is  a  column 
vector  whose  nth  component  is  the  value  of  the  incident  pressure  wave  at  a  reference 
point  on  Sn.  (For  acoustic  radiation  problems  Pinc  =  0.)  The  matrices  A  and  B  in 
equation  (4)  have  components  given  by 

Amn  =  \Smn~  !  dS{<j)  (5) 

Bmn  =  iurp  f  G(£m,(r)  dS{cr),  (6) 

Js„ 
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where  Smn  is  the  Kronecker  delta,  ( m  is  a  reference  point  on  Sm,  and  G((m,cr)  is  the 
free-space  Green’s  function 


G(Cm^) 

r{(m,v) 


e-ikr((.m,cr ) 
m>  l-r) 

the  distance  between  Cm  and  a. 


( k  —  u)/c ) 


The  norma!  derivative  of  G((m,cr)  is  given  by  the  expression 


dG((m,cr)  Gtfn,a)(.,  ,  1  ^r/_  ,  s  „,_m 

— a -  =  — 777 — 7  lk  +  77 - 7  K*  ~  Cm)  ■  n(a)\. 

On a  r2((m,cr )  \  r((,m,a)J 


(7) 


(8) 


It  is  well  known  that  the  solution  of  the  Helmholtz  integral  equation  for  p  in  terms 
of  v  and  p'nc  is  not  unique  at  certain  frequencies  corresponding  to  the  eigenfrequen- 
cies  of  the  interior  Dirichlet  problem.  In  the  neighborhood  of  these  eigenfrequencies 
the  matrix  A  in  equation  (4)  becomes  ill-conditioned.  To  overcome  this  numerical 
difficulty,  CHIEF  (Schenck,  1968)  supplements  the  equations  in  (4)  by  additional 
equations  of  the  same  form  in  which  each  is  a  point  interior  to  S.  The  resulting 
overdetermined  system  of  equations  has  a  unique  solution  in  the  least-squares  sense. 


The  total  acoustic  pressure  p  in  the  field  can  be  written  as  the  sum  of  the  incident 
pressure  p'nc  and  the  scattered  pressure  ps,  i.e., 

p(x)  =  ptnc{x)  +  pa{x).  (9) 

For  radiation  problems  p,nc  =  0  and  p  =  p’.  The  scattered  pressure  ps  at  a  field  point 
x  exterior  to  the  structure  can  be  approximated  by 

ps(x)  =  aT(x)P  +  bT  (x)V,  (10) 

where  a(x)  and  b(x)  are  column  vectors  having  components  given  by 


dS(a) 


/  x  [  dG(x,  a) 
an(x)  =  /  — ^ — —  aS(cr) 

Jsn  one 

bn(x )  =  iujp  /  G(x,  a)  dS(<r). 
Jsn 


(11) 

(12) 


In  the  farfield  of  S,  the  components  a„(x)  and  b„(x)  can  be  approximated  by 

'  /  —  t  kR. 

<■•(*)  =  (13) 

(14» 
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where  R  is  the  distance  from  an  origin  near  the  surface  S  to  the  field  point  x,  r (x) 
is  a  unit  vector  at  the  origin  in  the  direction  of  x,  and  6(cr)  is  the  position  vector  of 
the  integration  point  a  relative  to  the  origin  . 

Equations  (4)  and  (10)  can  be  combined  to  give 

Pa(x)  =  qT(x)V  +  Prs(*),  (15) 

where 

qT(x )  =  aT(x)A~1  B  +  bT(x) 

and 

pr3(x)  =  aT{x)A~1P,nc. 

It  can  be  seen  from  equation  (15)  that  prs(x)  is  the  rigid  scattering  from  the  body 
(i.e.,  the  scattered  pressure  when  V  is  constrained  to  be  zero). 

Applying  the  assumption  that  p  is  piecewise  constant  on  the  subdivisions  of  S  to 
equation  (2),  gives 

FTad  =  ~CDP,  (16) 

where 


and 

D  =  dia.g(Sx,S2,...,SN). 

Equations  (1)  and  (16)  can  be  combined  to  give 

U  =  ~{-u2M  +  KY'CDP  +  (— w2M  +  K)~xFd.  (17) 

Since  different  interpolation  schemes  are  used  in  the  finite  element  model  of  the 
structure  and  the  CHIEF  model  of  the  acoustic  field,  it  is  impossible  to  enforce  exact 
continuity  of  normal  displacement  across  S.  However,  this  continuity  is  approximately 
enforced  by  equating  the  CHIEF  normal  velocity  vn  to  the  average  of  the  finite  element 
normal  velocity  over  Sn,  i.e., 


v 


n 


iuu  ■  n  dS. 


Combination  of  equations  (3), (17),  and  (18)  gives 


V  =  iuCTU 

=  -iu.CT(-ui2M  +  KY'CDP  +  ivCri-u,1M  +  K)~'Fd. 


(18) 


(19) 
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Let  E  be  a  matrix  whose  columns  are  the  in  vacuo  normal  modes  of  the  structure, 

i.e., 

K  E  =  MEQ,  (20) 

where  Q  =  diag(u;j , . .  .  is  the  diagonal  matrix  of  eigenfrequencies.  Since  the 

modes  are  M-orthogonal,  they  can  be  normalized  so  that 


EtME  =  I.  (21) 

It  is  easily  verified  that  the  inverse  matrix  (—u>2M  +  A)-1  can  be  expressed  in  terms 
of  the  normal  mode  matrix  E  as  follows: 

(~u2M  +  A')'1  =  £(-u>2/  +  ft)-1  ET.  (22) 

This  form  is  convenient  since  the  diagonal  matrix  (— w2/  +  Q)  is  easily  inverted  at  all 
frequencies. 

Combination  of  equations  (19)  and  (22)  gives 

v  =  -iuCT  E{-u2 1  +  n)~1ETCDP  +  ilocte(-lj2i  +  nylETFD.  (23) 

Equations  (4)  and  (23)  can  now  be  combined  to  give 

HP  =  Ptnc  +  iu>BCT  E(—u2l  +  n)-'ETFD,  (24) 


where 

H  =  A  +  iu>BCT E(—u>2 1  +  \1)  xErCD. 

Equation  (24)  can  be  solved  for  the  surface  pressure  vector  P.  Once  P  is  determined, 
the  velocity  vector  V  can  be  obtained  from  equation  (23)  and  the  scattered  pressure 
p3(x)  can  be  computed  at  any  exterior  field  point  x  using  equation  (10)  or  equation 
(15). 
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FREQUENCY  INTERPOLATION 


In  this  section  certain  scaled  acoustic  quantities  will  be  constructed  which  vary 
slowly  with  frequency.  These  scaled  quantities  can  be  computed  with  a  coarser  fre¬ 
quency  resolution  than  is  desired  in  the  final  frequency  sweep,  and  the  intermediate 
frequency  values  can  then  be  obtained  by  interpolation.  The  only  acoustic  quantities 
involved  in  the  system  of  equations  (24)  for  P  are  the  matrices  A  and  B.  It  can 
be  seen  from  equations  (5)  to  (8)  that  the  components  of  A  and  B  both  involve  the 
quantity  e~,kT  which  varies  rapidly  with  frequency  when  the  separation  is  large.  This 
rapid  variation  in  the  integrands  of  Amn  and  Bmn  can  be  eliminated  by  factoring  out 
the  quantity  e~lkrrnn  where  rmn  is  the  distance  between  the  reference  points  on  Sm 
and  Sn.  Thus,  the  matrices  A  and  B  defined  by 

.in  =  e,krmn  f  dG^m'aldS{a)  (25) 

Jsn  one 

Bmn  =  iupe'kT™  [  G((m,cr )  dS(a )  (26) 

Js„ 

are  slowly  varying  with  frequency  and  are  simply  related  to  A  and  B  by 

A  —  LA  _  -ikrmn  A 

/*-mn  —  2Umn  c  rxmn  \^ '  ) 

Bmn  =  t-ikr™Bmn  .  (28) 

Linear  frequency  interpolation  can  be  used  on  A,  B  and  the  matrices  A,  B  can  be 
recovered  using  equations  (27)  and  (28). 


The  far-field  calculations  involve  the  acoustic  quantities  an(x),  bn(x)  defined  by 
equations  (13)  and  (14).  If  6(a)  in  these  equations  is  approximated  by  6n  —  £n  —  0, 
then  equations  (13)  and  (14)  reduce  to 


an(.f)  = 
bn  (  3-  ) 


r^Tr‘**M'fn  f  *(x)-nMdSM 

47t  R  Js„ 

iujpe  ikruv6nq 


(29) 

(30) 


For  equal  frequency  increments  the  quantity  e  ,krmn  occurring  in  equations  (27)  and 
(28)  and  the  quantity  occurring  in  equations  (29)  and  (30)  can  be  computed 

recursively  using 


ei(k+Ak)r  _  jkr  _  eiAkr 


(31) 


The  interpolation  scheme  can  now  be  summarized  as  follows: 
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1.  The  matrices  A  and  B  defined  in  equations  (25)  and  (26)  are  computed  at  a 
set  of  frequencies  which  have  a  wider  spacing  than  the  intended  final  frequency 
sweep. 

2.  The  integrals  fSn  r(x)  ■  n(a)  dS(a)  occurring  in  equation  (29)  are  computed  for 
each  subdivision  Sn  and  each  desired  far-field  point  x. 

3.  The  exponentials  e~,kTmn  and  e,k^^'^n  are  computed  for  a  starting  wavenum¬ 
ber  k0  and  a  wavenumber  increment  A k  corresponding  to  the  final  frequency 
resolution  desired. 

4.  At  each  wavenumber  k„  =  ko  +  (n  —  1)A k  the  matrices  A  and  B  are  computed 
using  linear  frequency  interpolation. 

5.  The  quantities  e~,kT,nn  and  e,k*(r)  are  updated  using  the  recursion  relation 
(31). 

6.  The  matrices  A  and  B  corresponding  to  kn  are  computed  using  equations  (27) 
and  (28). 

7.  Equation  (24)  is  solved  for  the  surface  pressures  P. 

8.  The  surface  normal  velocities  V  are  computed  using  equation  (23). 

9.  The  far-field  quantities  an(x)  and  bn(x)  are  computed  using  the  equations  (29) 
and  (30). 

10.  The  far-field  pressure  p3(x)  is  computed  for  each  field  point  x  using  equation 

(10). 

Examples  will  be  shown  in  the  next  section  which  illustrate  the  accuracy  and 
improvement  in  computational  efficiency  that  can  be  obtained  using  this  scheme. 
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NUMERICAL  RESULTS 


In  this  section  two  examples  will  be  presented  which  illustrate  the  accuracy 
and  improvement  in  computational  efficiency  that  can  be  obtained  using  frequency 
interpolation.  The  first  example  involves  the  backscattering  from  a  spherical  shell. 
The  ratio  of  the  thickness  to  the  mean  radius  of  the  shell  is  0.03  .  The  material 
parameters  of  the  shell  are 


•  Density  =  7669  Kg/m3 

•  Young’s  modulus  =  2.07  x  10n  N/m2 

•  Poisson’s  ratio  =  0.3  . 


The  material  parameters  of  the  water  are 


•  Density  =  998  Kg/m3 

•  Sound  Speed  =  1486  m/s. 


Figure  1  shows  the  backscattered  form  function  for  this  spherical  shell  versus  the 
normalized  frequency  ka  for  increasingly  coarse  interpolation  in  ka.  The  normalized 
frequency  ka  is  defined  by  ka  =  2irfa/c  where  /  is  the  frequency,  a  is  the  mean 
radius,  and  c  is  the  sound  speed  in  water.  The  ka  increment  Aka  used  in  the  acoustic 
calculations  (CHIEF)  was  taken  successively  to  be  0.01,  0.1,  0.2,  0.3.  It  can  be  seen 
from  figure  1  that  there  was  only  a  modest  decrease  in  accuracy  as  Aka  was  increased. 
Table  1  shows  the  number  of  frequencies  at  which  the  CHIEF  calculations  needed  to 
be  made  for  each  value  of  Aka. 

Table  1.  Number  of  CHIEF  frequencies 
required  for  sphere  backscattering  problem. 


Aka 

CHIEF  Frequencies 

0.01 

141 

0.10 

15 

0.20 

8 

0.30 

6 
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As  would  be  expected,  the  decrease  in  the  number  of  CHIEF  calculations  results 
in  a  large  reduction  in  overall  computation  time.  Figure  2  shows  the  central  pro¬ 
cessing  unit  (CPU)  times  required  on  a  CONVEX  CS-1  min '.supercomputer  for  100 
frequencies  and  91  far-field  angles  as  a  function  of  maximum  ka.  The  top  curve  shows 
the  CPU  time  required  for  the  acoustic  calculations  with  no  frequency  interpolation 
(Aka  =  .01).  The  next  lower  curve  gives  the  total  (acoustic,  finite-element,  and  cou¬ 
pling)  CPU  time  required  when  Aka  =  0.2  .  The  remaining  curves  show  the  CPU 
times  for  various  parts  of  the  calculation.  With  the  new  interpolation  scheme,  the 
time  is  dominated  at  low  ka  by  the  portion  of  the  program  which  couples  the  CHIEF 
and  finite  element  results.  At  higher  ka  the  finite  element  computations  dominate. 
In  general,  the  frequency  interpolation  scheme  has  resulted  in  an  overall  reduction  in 
computation  time  of  about  15:1. 
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Figure  1.  Backscattering  from  spherical  shell  utilizing  frequency  interpolation. 


Log  Convex  CPU  hours  for  100  frequencies 


Figure  2.  CPU  time  versus  ka  -  with  and  without  frequency  interpolation. 
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The  second  example  involves  the  backscattering  off  the  end  of  a  capped  cylindrical 
shell.  The  ratio  of  the  mean  radius  a  to  the  length  L  of  the  shell  is  0.125,  and  the 
ratio  of  the  thickness  h  to  the  mean  radius  a  of  the  shell  is  0.01016  .  The  ends  of  the 
cylinder  are  capped  with  circular  disks  having  the  same  thickness  as  the  shell.  The 
material  properties  of  both  the  shell  and  endcaps  are 

•  Density  =  8977  Kg/m3 

•  Young’s  modulus  =  2.09  x  1011  N/m2 

•  Poisson’s  ratio  =  0.308  . 


The  material  parameters  of  the  water  are 


•  Density  =  998  Kg/m3 

•  Sound  Speed  —  1486  m/s. 


Figure  3  shows  the  backscattered  normalized  pressure  for  an  incident  plane  wave 
travelling  along  the  axis  of  the  cylinder.  The  scattered  pressure  is  normalized  by  the 
high-frequency  plane-wave  approximation  to  the  rigid  backscattered  pressure.  The 
plane-wave  approximation  involves  setting  p,  =  pcv,  on  the  surface  where  p,  and 
v,  are  the  scattered  pressure  and  normal  velocity  respectively.  For  a  rigid  body  the 
scattered  and  incident  normal  velocities  on  the  surface  are  related  by  v,  =  —Vinc. 
Thus,  the  plane  wave  approximation  reduces  to  setting  p,  =  —pcvinc  on  the  surface. 
The  far-field  scattered  pressure  can  be  obtained  from  the  approximate  values  of  p„ 
and  v,  on  the  surface.  The  normalized  scattered  pressure  used  in  this  report  is  very 
similar  to  the  form  function  which  is  defined  to  be  the  magnitude  of  the  ratio  of  the 
scattered  pressure  to  the  geometric  acoustics  approximation  of  the  scattered  pressure. 
In  fact,  for  a  sphere,  the  two  normalizations  are  exactly  the  same.  The  normalized 
frequency  ka  is  defined  by  ka  =  2nfa/c  where  /  is  the  frequency,  a  is  the  mean 
radius,  and  c  is  the  sound  speed  in  water.  The  solid  curve  represents  computed 
values  at  ka  increments  of  0.025  with  no  interpolation.  The  dotted  curve  (which  is 
almost  indistinguishable  from  the  solid  curve)  corresponds  to  CHIEF  computations 
at  ka  increments  of  0.2  with  interpolation  to  0.025  increments. This  curve  represents 
CHIEF  computations  at  only  10  ka  values  out  of  a  total  of  71.  There  is  essentially  no 
loss  of  accuracy  in  this  frequency  range  for  this  degree  of  interpolation.  It  should  be 
noted  that  the  final  frequency  resolution  is  not  limited  to  the  value  ka  =  0.025  used 
in  this  example.  Interpolation  can  be  performed  to  any  desired  frequency  resolution 
without  any  further  CHIEF  computations. 
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Figure  4  shows  a  similar  comparison  for  CHIEF  ka  increments  of  0.4  .  This  rep¬ 
resents  only  six  CHIEF  frequency  computations  in  this  frequency  range.  Again  the 
agreement  is  very  good  showing  little  loss  of  accuracy  due  to  the  frequency  inter¬ 
polation.  Figure  5  shows  another  comparison  corresponding  to  ka  increments  of  1.0 
in  CHIEF.  Although  the  agreement  is  still  not  bad,  the  accuracy  has  deteriorated 
somewhat  in  the  upper  portion  of  the  frequency  band. 

Figure  6  illustrates  the  importance  of  scaling  the  components  of  the  CHIEF  A  and 
B  matrices  by  e-,fcrmn.  The  solid  curve  again  represents  calculations  at  ka  increments 
of  0.025  with  no  interpolation.  The  dotted  curve  represents  ka  increments  of  0.4 
for  the  CHIEF  computations  with  interpolation  performed  on  the  A  and  B  matrices 
directly  without  any  scaling.  It  is  clear  that  the  agreement  is  not  nearly  as  good  as 
it  was  for  the  scaled  interpolation  case  shown  in  figure  4. 
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Figure  3.  Backscattering  from  capped  cylinder,  CHIEF  ka  increment  of  0.2. 


m 


Figure  4.  Backscattering  from  capped  cylinder,  CHIEF  ka  increment  of  0.4. 
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Figure  5.  Backscattering  from  capped  cylinder,  CHIEF  ka  increment  of  1.0. 


Figure  6.  Effect  of  scaling  on  frequency  interpolation. 
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CONCLUSIONS  AND  RECOMMENDATIONS 


For  coupled  structural  acoustic  scattering  problems,  the  spectral  response  often 
varies  rapidly  over  a  large  dynamic  range,  requiring  the  response  to  be  calculated  at 
many  closely  spaced  frequencies.  In  such  problems,  the  computation  time  is  often 
dominated  by  the  calculation  of  acoustic  interactions  which  do  not  vary  rapidly  with 
frequency.  It  has  been  shown  in  this  report  that  large  savings  in  computation  time 
can  be  achieved  without  significant  loss  of  accuracy  by  interpolating  on  appropriately 
scaled  intermediate  acoustic  quantities.  These  quantities  are  computed  on  a  coarse 
frequency  grid  and  then  interpolated  to  any  desired  frequency  resolution  in  order  to 
obtain  the  final  response.  The  overall  savings  in  computation  time  is  typically  more 
than  an  order  of  magnitude. 

It  is  recommended  that  this  frequency  interpolation  technique  be  incorporated  in 
other  similar  computer  codes  dealing  with  structural  acoustic  scattering  problems. 
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